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Introduction 

• Dual numbers extend real numbers, 
similar to complex numbers. 

• Complex numbers adjoin an element /, for 
which i 2 = -1. 

• Dual numbers adjoin an element e, for 
which e 2 = 0. 
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Complex Numbers 

• Complex numbers have the form 

z = a + b i 

where a and b are real numbers. 

• a = real(z) is the real part, and 

• b = imag(z) is the imaginary part. 
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Complex Numbers (cont'd) 

• Complex operations pretty much follow 
rules for real operators: 

• Addition: 

(a + b i) + (c + d /) = (a + c) + (b + d) i 

• Subtraction: 

(a + b /) - (c + d /) = (a - c) + (b - d) i 
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Complex Numbers (cont'd) 

• Multiplication: 

(a + b /) (c + d /) = (ac - bd ) + ( ad + be) i 

• Products of imaginary parts feed back into 
real parts. 
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Dual Numbers 

• Dual numbers have the form 

z = a + b £ 

similar to complex numbers. 

• a = real(z) is the real part, and 

• b = dual(z) is the dual part. 
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Dual Numbers (cont'd) 

• Operations are similar to complex 
numbers, however since e 2 = 0, we have: 

(a + b e) (c + d e) = (ac + 0) + (ad + bc)£ 

• Dual parts do not feed back into real 
parts! 




GAME DEVELOPERS CONFERENCE' 2013 



MARCH 25-29, 2013 GDCONF.COM 



Dual Numbers (cont'd) 

• The real part of a dual calculation is 
independent of the dual parts of the inputs. 

• The dual part of a multiplication is a 
"cross" product of real and dual parts. 
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Taylor Series 

• Any value f(a + h) of a smooth function f 
can be expressed as an infinite sum: 

. f\d) 7 f\a) , 2 

f ( CL + //) — f ( Cl ) H h H h + * * * 

1 ! 2 ! 

where f' f f" r A n ) are the first, second, ..., 
/7-th derivative of f. 
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Taylor Series Example 
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Taylor Series and Dual Numbers 

• For f(a + b e), the Taylor series is: 

fia + bs) = f(ay + nfbs + ... O 

• All second- and higher-order terms 
vanish! 

• We have a closed-form expression that 
holds the function and its derivative. 
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Real Functions on Dual Numbers 

• Any differentiable real function f can be 
extended to dual numbers, as: 

f(a + b £) = f(a ) + b f\a ) s 

• For example, 

sin{a + b s) = sin{a ) + b cos(a) e 
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Automatic Differentiation 

• Add a unit dual part to the input value of 
a real function. 

• Evaluate function using dual arithmetic. 

• The output has the function value as real 
part and the derivate's value as dual part: 



f{a + s) = f(a ) + f\a) e 
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How does it work? 

• Check out the product rule of 
differentiation: (/ -gy=f-g'+f'-g 

• Notice the "cross" product of functions and 
their derivatives. 

• Recall that 

(a + a'e){b + b'e) = ab + ( ab'+ a'b)e 
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Automatic Differentiation in C+ + 

• We need some easy way of extending 
functions on floating-point types to dual 
numbers... 

• ...and we need a type that holds dual 
numbers and offers operators for 
performing dual arithmetic. 
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Extension by Abstraction 

• C++ allows you to abstract from the 
numerical type through: 

• Typedefs 

• Function templates 

• Constructors and conversion operators 

• Overloading 

• Traits class templates 
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Abstract Scalar Type 

• Never use built-in floating-point types, 
such as float or double, explicitly. 

• Instead use a type name, e.g. scalar, 
either as template parameter or as typedef, 



typedef float Scalar; 
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Constructors 

• Built-in types have constructors as well: 

• Default: float () == 0 . Of 

• Conversion: float (2) == 2. Of 

• Use constructors for defining constants, 
e.g. use scalar (2) , rather than 2 . Of or 

(Scalar) 2 . 
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Overloading 

• Operators and functions on built-in types 
can be overloaded in numerical classes, 

SUCh as std: : complex. 

• Built-in types support operators: + 

• ...and functions: sqrt, pow, sin, ... 

. NB: Use <cmath> rather than <math.h>. 
That is, use sqrt NOT sqrtf on floats. 
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Traits Class Templates 

• Type-dependent constants, such as the machine 
epsilon, are obtained through a traits class defined 

in <limits > . 

• Use std : : numeric_limits<Scalar> : : epsilon ( ) 

rather than flt epsilon in C+ + . 

• Either specialize std: :numeric_iimits for your 
numerical classes or write your own traits class. 
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Example Code (before) 

float smoothstep ( float x) 

{ 

if (x < 0 . Of ) 
x = 0 . 0 f ; 

else if (x > 1 . Of ) 
x = 1 . 0 f ; 

return (3 . Of - 2 . Of * x) * x * x; 

} 




GAME DEVELOPERS CONFERENCE' 2013 



MARCH 25-29, 2013 GDCONF.COM 



Example Code (after) 

template <typename T> 

T smoothstep (T x) 

{ 

if (x < T() ) 
x = TO; 

else if (x > T(l)) 
x = T ( 1 ) ; 

return (T(3) - T(2) * x) * x * x 

} 



> • 
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Dual Numbers in C+ + 

• C++ has a standard class template 
std : : compiex<T> for complex numbers. 

• We create a similar class template 
Duai<T> for dual numbers. 

• Duai<T> defines constructors, accessors, 
operators, and standard math functions. 
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Dual<T> 

template <typename T> 
class Dual 
{ 

• • • 

private : 

T mReal; 

T mDual; 

}; 
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Dual<T>: Constructor 

template <typename T> 

Dual<T> : : Dual (T real = T(), T 
: mReal (real) 
f mDual (dual ) 

{ } 



dual = T()) 



Dual<Scalar> zl; // zero initialized 
Dual<Scalar> z2(2); // zero dual part 

Dual<Scalar> z3(2, 1); 
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Dual<T>: operators 

template <typename T> 

Dual<T> operator* (Dual<T> a, Dual<T> b) 

{ 

return Dual<T> ( 

a. real () * b.real(), 

a. real () * b.dual() + 

a. dual () * b.real() 

) ; 



} 
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Dual<T>: Standard Math 

template <typename T> 

Dual<T> sqrt (Dual<T> z) 

{ 

T tmp = sqrt ( z . real ()) ; 
return Dual<T> ( 

tmp, 

z.dual() * T(0.5) / tmp 

) ; 



} 
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Curve Tangent 

• For a 3D curve 

pO) = 0(0, yiO, z(f )), where t e [a, b\ 

The tangent is 

P'OO 

IIpTOII ’ 



where p'(0 = (*'(0, y'(t ) , z'(0) 
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Curve Tangent 



• Curve tangents are often computed by 
approximation: 



pOi)-pOo) 

||pOi)-p(* 0 ) 



9 



where f, = t 0 + h 



for tiny values of h. 
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Curve Tangent 

p(t 0 ) , £ 
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Actual 

tangent 



p(ti) 
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Curve Tangent: Bad #2 



t x drops outside 
parameter domain 
(ti > b) 
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Curve Tangent: Duals 

• Make a curve function template using a 
class template for 3D vectors: 

template <typename T> 
Vector3<T> curveFunc (T x) 



• < 
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Curve Tangent: Duals (cont'd) 

• Call the curve function using a dual 

number X = Dual<Scalar> (t, 1), 

(add £ to parameter t) : 



Vector3<Dual<Scalar> > y = 

curveFunc (Dual<Scalar> (t, 1) ) ; 
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Curve Tangent: Duals (cont'd) 

• The real part is the evaluated position: 

Vector3<Scalar> position = real (y) 

• The normalized dual part is the tangent at 
this position: 

Vector3<Scalar> tangent = 
normalize (dual (y) ) ; 



• < 
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Line Geometry 

• The line through points p and q can be 
expressed explicitly as: 

x(t) = P + (q - P )t, and 

• Implicitly, as a set of points x for which: 

(q - p) x x + p x q = 0 
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Line Geometry 




p x q is orthogonal to the plane opq, and its length is equal 
to the area of the parallellogram spanned by p and q 
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Line Geometry 




All points x on the line pq span with q - p a parallellogram 
that has the same area and orientation as the one spanned 
by p and q. 
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Plucker Coordinates 

• Plucker coordinates are 6-tuples of the 
form (u x , u y , u z/ v x/ v y/ v z ), where 

u = (u x , u yf u z ) = q - p, and 

v = (v x , v y , v z ) = p x q 
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Plucker Coordinates (cont'd) 

• For (u^vj and (u 2 :v 2 ) directed lines, if 

u 1 • v 2 + y 1 • u 2 is 

zero: the lines intersect 

positive: the lines cross right-handed 
negative: the lines cross left-handed 
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Triangle vs. Ray 




If the signs of permuted dot products of the ray 
and edges are all equal, then the ray intersects the 
triangle. 
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Plucker Coordinates and Duals 

• Dual 3D vectors conveniently represent 
Plucker coordinates: 

Vector3<Dual<Scalar> > 



• For a line (u:v), u is the real part and v is 
the dual part. 
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Dot Product of Dual Vectors 

• The dot product of dual vectors u ± + v x £ 
and u 2 + v 2 £ is a dual number z, for which 

real(z) = u x • u 2/ and 
dual(z) = u x • v 2 + v x • u 2 

• The dual part is the permuted dot product 
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Angle of Dual Vectors 



• For a and b dual vectors, we have 



6 + ds = arccos 



( 


a *b ^ 


V 


a 


b U 



where 9 is the angle and d is the signed 
distance between the lines a and b. 
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Translation 

• Translation of lines only affects the dual 
part. Translation of line pq over c gives: 

• Real: (q + c) - (p + c) = q - p 

• Dual: (p + c) x (q + c) 

= p x q + c x (q - p) 

• q - p pops up in the dual part! 
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Rotation 

• Real and dual parts are rotated in the 
same way. For a rotation matrix R: 

• Real: Rq - Rp = R(q - p) 

• Dual: Rp x Rq = R(p x q) 

• The latter holds for rotations only! That is, 
R performs no scaling or reflection. 
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Rigid-Body Transform 



• For rotation matrix R and translation vector c, 
the dual 3x3 matrix M with 
real(M) = R, and 



dual(M) = [c]„R 




o -c 



X 



C 



o 



R 



maps Plucker coordinates to the new reference 
frame. 




GAME DEVELOPERS CONFERENCE' 2013 



MARCH 25-29, 2013 GDCONF.COM 



Screw Theory 

• A screw motion is a rotation about a line 
and a translation along the same line. 

• " Any rigid body displacement can be 
defined by a screw motion (Chasles) 




GAME DEVELOPERS CONFERENCE' 2013 



MARCH 25-29, 2013 GDCONF.COM 



Chasles' Theorem (Sketchy Proof) 

• Decompose translation into a term along 
the line and a term orthogonal to the line. 

• Translation orthogonal to the axis of 
rotation offsets the axis. 

• Translation along the axis does not care 
about the position of the axis. 
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Translations Orthogonal to Axis 
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Example: Rolling 
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Ball 
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Dual Quaternions 

• Unit dual quaternions represent screw 
motions. 

• The rigid body transform over a unit 
quaternion q and vector t is: 



q h — tq<? 
2 



Here, t is a quaternion 
with zero scalar part. 
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Where is the Screw? 



• A unit dual quaternion can be written as 



cos 



r 



v 



0 + ds 



\ 



J 



sin 



r 



v 



0 + ds 



\ 



J 



(u + v^) 



where 9 is the rotation angle, d, the 
translation distance, and u + vs, the line 
given in Plucker coordinates. 




GAME DEVELOPERS CONFERENCE' 2013 



MARCH 25-29, 2013 GDCONF.COM 



Rigid-Body Transform Revisited 

• Similar to 3D vectors, Plucker coordinates 
can be transformed using dual quaternions. 

• The mapping of a dual vector v according 
to a screw motion q is 



v' = qvq* 
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Traditional Skinning 

• Bones are defined by transformation 
matrices T, relative to the rest pose. 

• Each vertex is transformed as 



p - A T iP + • • • + Wp - + • • • + \T n )p 

Here, A, are blend weights. 
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Traditional Skinning (cont'd) 

• A weighted sum of matrices is not 
necessarily a rigid-body transformation. 

• Most notable artifact is "candy wrapper": 
The skin collapses while transiting from one 
bone to the other. 




GAME DEVELOPERS CONFERENCE' 2013 



MARCH 25-29, 2013 GDCONF.COM 



Candy Wrapper 
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Dual Quaternion Skinning 

• Use a blend operation that always returns 
a rigid-body transformation. 

• Several options exists. The simplest one 
is a normalized lerp of dual quaternions: 

_ AQi + • • • + 

X l q l + • • • + A n q n 
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Dual Quaternion Skinning (cont'd) 

• Can the weighted sum of dual quaternions 
ever get zero? 

• Not if all dual quaternions lie in the same 
hemisphere. 

• Observe that q and -q are the same 
pose. If necessary, negate each q, to dot 
positively with q 0 . 
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Further Uses 

• Motor Algebra: Linear and angular 
velocity of a rigid body combined in a dual 
3D vector. 

• Spatial Vector Algebra: Featherstone 
uses 6D vectors for representing velocities 
and forces in robot dynamics. 
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Conclusions 

• Abstract from numerical types in your 
C++ code. 

• Differentiation is easy, fast, and exact 
with dual numbers. 

• Dual numbers have other uses as well. 
Explore yourself! 
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Thank You! 

• For sample code, check out free* MoTo 
C++ template library on: 

https : / / code . google . com/p /mot ion- tool kit/ 



(*) gratis (as in "free beer") and libre (as in "free speech") 




